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Abstract: We present the first determination of the strong coupling constant from a fit of 
next-to-next-to-leading order QCD predictions to event-shape variables, measured in e + e~ 
annihilations at LEP. The data have been collected by the ALEPH detector at centre-of- 
mass energies between 91 and 206 GeV. Compared to results of next-to- leading order fits 
we observe that the central fit values are lower by about 10%, with considerably reduced 
scatter among the results obtained with different event-shape variables. The dominant 
systematic uncertainty from renormalization scale variations is reduced by a factor of two. 
By combining the results for several event-shape variables and centre-of-mass energies, we 
find 

a 8 (Mf) = 0.1240 ± 0.0008 (stat) ± 0.0010 (exp) ± 0.0011 (had) ± 0.0029 (theo). 

Keywords: QCD, Jets, LEP Physics, NLO and NNLO Computations, strong coupling 
constant. 
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1. Introduction 

Quantum Chromodynamics (QCD) is generally accepted to be the correct theory for the 
description of strong interactions between quarks and gluons [1]. If the quark masses and 
mixing angles are fixed, then the only free parameter of the theory is the strong coupling 
constant, a s . Therefore it is of paramount importance to measure this parameter at the best 
possible precision. In particular, measurements based on different underlying processes, at 
different energy scales, constitute an important consistency check of the theory and are used 
to prove the so-called running of the coupling constant, ie., the decrease of the coupling 
strength with increasing energy scale [2] . During the last twenty years an enormous wealth 
of measurements has become available, using many different processes with different initial 
and final states [3,4]. These measurements resulted in values for a s (Q 2 ), at different 
scales Q, which are perfectly consistent with the expected running of the coupling. When 
evolved to a single scale, most commonly chosen to be the mass of the Z boson, Mz, the 
measurements agree with each other within a few percent, and thus can be used to extract 
a world-average value. The particle data group [4] quotes a s (M|) = 0.1176±0.0020, which 
has been improved upon in a more recent global analysis to a s (M|) = 0.1189 ± 0.0010 [5]. 

Several classes of observables can be identified as useful for an a s measurement. In- 
clusive observables, such as the ratio of the cross sections for e + e~ — > hadrons and 
e + e~ — > or sum rules in deep-inelastic scattering, do not resolve the structure of 
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the hadronic final state. The theoretical predictions are obtained with perturbative QCD 
(pQCD) and typically known up to next-to-next-to-leading order (NNLO) in the strong 
coupling constant. Non-perturbative effects related to the transition (hadronisation) of 
partons (quarks and gluons) to hadrons are strongly suppressed by the scale Q of the 
process. Exclusive observables, such as event-shape variables in e + e _ annihilations, dis- 
cussed in detail below, resolve topological properties of the hadronic final state and thus 
are sensitive to gluon radiation. This class of observables shows stronger sensitivity to 
hadronisation effects, at least in phase-space regions characterised by soft and collinear 
gluon radiation. Finally, spectroscopic properties of mesons and baryons, such as masses or 
energy level splittings of heavy quark-antiquark resonances, are necessarily determined by 
non-perturbative phenomena. Major recent progress in calculations within non-relativistic 
QCD and lattice QCD has resulted in very precise a s determinations from this last type 
of observables [6]. 

Thanks to the clean initial state and the very large event statistics, e + e~ annihila- 
tions to hadrons at the LEP collider at centre-of-mass energies between 91 and 206 GeV 
constituted an excellent laboratory for precision tests of QCD in general, and precise mea- 
surements of the strong coupling constant in particular [7-10]. These are complemented by 
data obtained by the SLD experiment [11] at SLAC at centre-of-mass energy of 91 GeV and 
by data obtained at the DESY PETRA collider at lower energies, especially the reanalysis 
of data from the JADE experiment [12]. 

Many of those measurements were based on event-shape variables such as thrust [13]. 
These variables probe the structure of the hadronic final state. In the leading-order picture, 
e + e~ annihilation to hadrons occurs via e + e~ — > qq and subsequent hadronisation to stable 
hadrons, resulting in a back-to-back two-jet structure of the event. At the next order 
of perturbation theory, gluon radiation off quarks will lead to deviations from this two- 
jet structure due to the appearance of additional jets. In general, event-shape variables 
quantify the structure of an event by a single measure. 

Until recently, pQCD predictions for event-shape observables were known up to next- 
to-leading order (NLO) in the strong coupling. As a consequence, the residual sensitivity to 
the choice of the unphysical renormalisation scale, fi, resulted in systematic uncertainties 
of the a s measurements at the 10% level. Those systematic uncertainties are typically 
estimated by varying jj, over a broad range, eg., 1/2 < fi/Q < 2, since this variation probes 
most of the missing higher-order contributions. Using a very large set of event shapes, the 
DELPHI collaboration observed a large spread in the a s (M|) values extracted from the 
different observables [14], when fixing fi = Mz- In addition, the measurements preferred 
rather large values of a s (M|) ~ 0.125 — 0.13 compared to the world average. On the other 
hand, simultaneously fitting a s (M|) and the scale \x resulted in a considerably reduced 
spread of a s (M|) values, however, at the price of a large spread of preferred ji parameters, 
which in addition turned out to be rather small compared to the "natural" choice fi = Mz- 
Similar observations were made by the other LEP collaborations and on the SLD data [15]. 

An important part of the missing higher-order contributions was identified to be due to 
large logarithms in the event-shape variable, related to the incomplete cancellation of real 
and virtual corrections due to gluon radiation [16-21]. Once these logarithms are resummed 
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to all orders in pQCD and consistently combined with the fixed-order calculations, the re- 
sulting improved predictions lead to a better description of the event-shape variables, most 
notably close to the two-jet region of the phase space. Consequently, the a s measurements 
from fits of these improved predictions were characterised by reduced systematic uncer- 
tainties (now at the 5% level) and showed better agreement with other determinations, 
with typical values in the range a s (M|) = 0.118 — 0.125. Other systematic uncertainties, 
from experimental effects and most notably from the estimation of non-perturbative con- 
tributions, were smaller than those from the //-scale variations. Non-perturbative hadroni- 
sation corrections are estimated by either using phenomenological models (string or cluster 
fragmentation), implemented in Monte Carlo simulations [22-24], or by adding power-law 
corrections to the purely perturbative predictions [25]. 

From all this it became clear that major progress in the pQCD calculations was nec- 
essary in order to push the uncertainties of a s measurements from event shapes below the 
5% range, approaching the precision obtained with inclusive observables or lattice calcu- 
lations. Such major progress has recently been made with the calculation of the NNLO 
corrections to event-shape variables [26,27]. Consequently, it is of great interest to include 
these new corrections in fits for a s (M|). In this paper we describe first fits of this kind 
and discuss the results. The measurements are obtained with event-shape distributions 
measured by the ALEPH experiment [28] at centre-of-mass energies between 91 and 206 
GeV. We concentrate on the comparison of the results found at NNLO with those at NLO 
and NLO matched to resummation of leading and next-to- leading logarithms (NLLA) , as 
used in [28]. In a forthcoming publication results using NNLO matched to NLLA will be 
discussed in detail. 



2. Theoretical predictions 



The perturbative expansion for the distribution of a generic observable y up to NNLO at 
e + e~ centre-of-mass energy y/s, for a renormalisation scale fj 2 = s and a s = a s (s) is given 
by 



1 da 



dcr / 1 



a,\ dA 



Chad dy V2vr/ dy \2ttJ dy \2irJ dy 



(2.1) 



Here the event-shape distribution is normalised to the total hadronic cross section Chad- 
With the assumption of massless quarks we have at NNLO 
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where C3 = 1.202056903 . . ., and the QCD colour factors are 

N 2 - 1 1 
C A = N, C F = ^ r , T R = - (2.5) 

for N = 3 colours and N F light quark flavours. 

In practice, we use the perturbative coefficients A, B and C, computed in the MS- 
scheme in [27], which are all normalised to ctq: 

1 da fa s \ dA /a s \ 2 dB /a s \3 dC m , A . , 

t v = y n + (s) * + y * + o( "- ) • <2 - 6) 

However, A, and C are straightforwardly related to A, B and (7, 

A = A , 

B = B-^C F A, 

3„ „ /9^ 2 



C = C--C F B+ y-C F -K 2 j A. (2.7) 

These coefficients are computed at a renormalisation scale fixed to the centre-of-mass en- 
ergy, and depend therefore only on the value of the observable y. They explicitly include 
only QCD corrections with non-singlet quark couplings and are therefore independent of 
electroweak couplings. At 0{a 2 ), these amount to the full corrections, while the O(af) cor- 
rections also receive a pure-singlet contribution. This pure-singlet contribution arises from 
the interference of diagrams where the external gauge boson couples to different quark lines. 
In four-jet observables at 0{aP s ), these singlet contributions were found to be extremely 
small [29]. 

The QCD coupling constant evolves according to the renormalisation group equation, 
which is to NNLO: 



2 da s (/T 



(2.8) 



with the MS-scheme coefficients 
11 C A -4T R N F 



A = 



6 

17 C\ - 10 C A T R N F - 6 C F T R N F 
6 



h = — (2857C1 + 108 C 2 F T R N F - 1230 C F C A T R N F - 2830 C 2 A T R N F 

+264 CpT^Np + 316 C A T R N 2 F ) . (2.9) 

The above equation is solved by introducing A^j^f^ as integration constant with L = 
2 log(/z/A^?^), yielding the running coupling constant: 



a 
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In terms of the running coupling a s (fi 2 ), the NNLO expression for event-shape distri- 
butions becomes 



Chad dy 



1 da 




) 



(2.11) 



The coefficients A, B and C have been computed for several event-shape variables [26, 
27]. The calculation is carried out using a newly developed parton-level event generator 
programme EERAD3 which contains the relevant matrix elements with up to five external 
partons [30-33]. Besides explicit infrared divergences from the loop integrals, the four- 
parton and five-parton contributions yield infrared divergent contributions if one or two of 
the final state partons become collinear or soft. In order to extract these infrared diver- 
gences and combine them with the virtual corrections, the antenna subtraction method [34] 
was extended to NNLO level [35] and implemented for e + e~ — > 3jets and related event- 
shape variables [36] . The analytical cancellation of all infrared divergences serves as a very 
strong check on the implementation. EERAD3 yields the perturbative A, B and C coef- 
ficients as histograms for all infrared-safe event-shape variables related to three-particle 
final states at leading order. As a cross check, the A and B coefficients have also been 
obtained from an independent integration [37] of the NLO matrix elements [38], showing 
excellent agreement. While A and B can be computed to very high numerical accuracy, 
the computation of C is so CPU-intensive that we could obtain only results with numerical 
integration errors of typically a few percent. These numerical errors on the C-coefficient 
will be included in the analysis of a s as described below. 

For small values of y, the fixed-order expansion, eq. (2.11), fails to converge, because the 
fixed-order coefficients are enhanced by powers of ln(l/y), yA(y) ~ lny, yB(y) ~ ln 3 y and 
yC(y) ~ ln 5 y. In order to obtain reliable predictions in the region of y <C 1 it is necessary 
to resum entire sets of logarithmic terms at all orders in a s . A detailed description of the 
predictions at next-to- leading- logarithmic approximation (NLL A) can be found in Ref . [39] . 
There the matching of the NLO and NLLA expressions, which is necessary in order to avoid 
double counting, is also described. The extension of this matching procedure to NNLO 
with NNLA is outlined in [17] and will be discussed in detail in a separate publication. 

3. Observables and data sets 

We have studied the six event-shape distributions thrust T [13], heavy jet mass Mh [41], 
total and wide jet broadening (Bt,Bw) [42], C-parameter C [43] and the Durham jet 
resolution parameter y^ [44]. The definitions of these variables and a discussion of their 
properties can be found in Refs. [28] and [39,40]. 
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The measurements have been carried out by the ALEPH collaboration [28] 1 , at centre- 
of-mass energies of 91.2, 133, 161, 172, 183, 189, 200 and 206 GeV. Earlier measurements 
and complementary data sets from the LEP experiments and from SLD can be found 
in Refs. [7-11]. The event-shape distributions were computed using the reconstructed 
momenta and energies of charged and neutral particles. The measurements have been 
corrected for detector effects, ie., the final distributions correspond to the so-called particle 
(or hadron level). The particle level is defined by stable hadrons with a lifetime longer 
than 10 -9 s after hadronisation and leptons according to the definition given in [47]. In 
addition, at LEP2 energies above the Z peak they were corrected for initial-state radiation 
effects. At energies above 133 GeV backgrounds from WW, ZZ and Z7* were subtracted 
following the procedure given in [28]. 

Backgrounds, mainly from W-pair production, were subtracted. The experimental 
uncertainties were estimated by varying event and particle selection cuts. They are below 
1% at LEP1 and slightly larger at LEP2. For further details we refer to Ref. [28]. 

4. Determination of the strong coupling constant 

The coupling constant a s is determined from a fit of the perturbative QCD predictions to 
measured event-shape distributions. The procedure adopted here follows closely the one 
described in Ref. [28]. The six event-shape variables used are T, Mh, Bt, B\y, C and 
— lm/3. The perturbative predictions for the distributions, as described in section 2, are 
calculated to the same order of perturbation theory for all of these variables. The data set 
collected by ALEPH at the Z peak with high statistics allows for quantitative comparisons 
of different pQCD predictions and resulting fits for a s . The size of missing higher orders, 
which are inherently difficult to assess, can be different for different variables. Therefore, 
a combination of measurements using several variables yields a better estimator of a s 
than using a single variable. Furthermore, the spread of values of a s is an independent 
estimation of the theoretical uncertainty. At centre-of-mass energies above the Z peak the 
statistical uncertainties are larger and background conditions are more difficult compared 
to the peak data. Therefore a combination of measurements is particularly important for 
those energies. We adopt the same combination procedure as described in [28], which is 
based on weighted averages and takes into account correlations between the event-shape 
variables. At energies above Mz the measured statistical uncertainties are replaced by 
expected statistical errors, obtained from a large number of simulated experiments. This 
avoids biases towards lower values of a s in the case of small event statistics. The bias 
originates from larger weights of downward fluctuating bins in the distributions compared 
to upward fluctuations, as outlined in Ref. [28]. 

Event-shape distributions are fit in a central region of three-jet production, where a 
good perturbative description is available. The fit range is placed inside the region where 
hadronisation and detector corrections are below 25% and the signal-to-background ratio at 
LEP2 is above one. At the higher LEP2 energies the good perturbative description extends 

lr The tables with numbers and uncertainties for all variables can be found at 
http: //aleph. web. cern. ch/aleph/QCD/alephqcd. html. 
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further into the two-jet region, while in the four-jet region the background becomes large. 
Thus the fit range is selected as a result of an iterative procedure balancing theoretical, 
experimental and statistical uncertainties. The data are corrected for detector effects, 
for background from four-fermion processes and for the residual ISR contribution. The 
background from WW events increases with energy, and after subtraction the content of 
some bins of the distribution becomes negative. For this reason the fit range is restricted 
to a region with good signal-to-background ratio. 

Since very recently distributions of infrared- and collinear-safe observables at the par- 
ton level can be computed in perturbative QCD to third order in a s , as described in section 
2 and in more detail in Ref. [36]. Here we concentrate on fits of NNLO predictions and 
compare them to pure NLO and matched NLO+NLLA predictions as used in the analysis 
of Ref. [28] . The nominal value for the renormalisation scale = jjl/Q is unity. All these 
calculations neglect quark masses. Quark mass effects [45] are relevant for the b quark 
at \/s = Mz, where the fraction of bb events is large, while Q is still moderate. NLO 
calculations including a quark mass indicate that the expected change in a s is of the order 
of 1% at M z [46]. The effect is scaling with M 2 /Q 2 and decreases to 0.2-0.3% at 200 GeV. 
Mass corrections were computed to second order using the matrix elements of [48] . A pole 
b-quark mass Mb = 5 GeV/c 2 was used and Standard Model values were taken for the 
fraction of bb events. It is worth noting that, since these mass effects are known at NLO 
only, the NNLO predictions can be corrected only in a partial and not fully consistent 
manner. 

The perturbative QCD prediction is corrected for hadronisation and resonance decays 
by means of a transition matrix, which is computed with the Monte Carlo generators 
PYTHIA, HERWIG and ARIADNE, all tuned to global hadronic observables at M z [47]. 
The parton level is defined by the quarks and gluons present at the end of the parton shower 
in PYTHIA and HERWIG and the partons resulting from the color dipole radiation in 
ARIADNE. Albeit the generators' limitation to leading-order QCD, a reasonable agreement 
is observed between the generator parton level and the NNLO prediction inside the fit 
range within a few percent, as shown in Fig. 1. Corrected measurements of event-shape 
distributions are compared to the theoretical calculation at particle level. 

The value of a s is determined at each energy using a binned least-squares fit. The fit 
programs of Ref. [28] have been extended to incorporate the NNLO calculations. It has 
been verified that at NLO+NLLA and using the fit range advocated by ALEPH the results 
on a s reported in Ref. [28] are reproduced. Only statistical uncertainties arising from the 
limited number of observed events, from the number of simulated events used to calculate 
hadronisation and detector corrections and from the integration procedure of the NNLO 
coefficient functions are included in the x 2 of the fit. Its quality is good for all variables at 
all energies. Nominal results, based on (2.1) and using the fitted values of a s are shown 
in Figs. 2 and 3 together with the measured distributions. The resulting measurements of 
a s (Q) are given in Table 1 for 91.2 to 172 GeV and in Table 2 for 183 to 206 GeV. 

Furthermore, in Figs. 4 and 5 and in Table 3 we compare the results found when 
fitting different perturbative predictions to the high-precision data at Mz- In particular 
we compare the NNLO fits to those obtained with NLO and NLO+NLLA predictions. A 
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Figure 1: Comparison of the parton- level distributions predicted by the generators and the NNLO 
calculation for the heavy jet mass (left) and -ln(y3) (right). 

detailed description of the latter is given in Ref. [28]. As shown in Figs. 4 and 5 a good 
description by the NNLO predictions is found for an extended range over the perturbative 
region compared to NLO. However, in the two-jet region the resummed prediction still 
yields a better result. Therefore, compared to Ref. [28], the fit range for NNLO and NLO 
is slightly shifted into the perturbative region. We observe a clear improvement in the fit 
quality when going to NNLO accuracy. Compared to NLO the value of a s is lowered by 
about 10%, but still higher than for NLO+NLLA. Systematic theoretical and experimental 
uncertainties are discussed in the following section. 

5. Systematic Uncertainties of a. s 
5.1 Experimental Uncertainties 

A detailed description of the experimental systematic uncertainties on the measured event- 
shape distributions is found in Ref. [28]. The resulting uncertainties on the fitted a s values 
are estimated in a way similar to that for the event shapes themselves. Changes of the 
distributions under variations of event and particle selection cuts lead in general to small 
changes in a s . In the fit procedure the same expected statistical uncertainties (cf. section 
4) are assumed for all variants of the distribution. This procedure reduces purely statistical 
components in the systematic effect, which are potentially large at LEP2 energies. The 
total experimental systematic uncertainties of a s at LEP2 are between 0.5% and 1.5%, 
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Figure 2: Distributions measured by ALEPH, after correction for backgrounds and detector 
effects, of thrust, heavy jet mass, total and wide jet broadening at energies between 91.2 and 206 
GcV together with the fitted NNLO QCD predictions. The error bars correspond to statistical 
uncertainties. The fit ranges cover the central regions indicated by the solid curves, the theoretical 
predictions extrapolate well outside the fit ranges, as shown by the dotted curves. The plotted 
distributions are scaled by arbitrary factors for presentation. 



dominated by limitations of the correction scheme for initial state radiation. Those at 
LEP1 are below 1% and dominated by imperfections of the simulation of reconstructed 
neutral hadronic energy deposits. As expected, the experimental uncertainties determined 
when fitting the NNLO prediction are generally the same as reported by ALEPH using 
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Figure 3: Distributions measured by ALEPH, after correction for backgrounds and detector ef- 
fects, of C-parameter and the three-jet transition variable at energies between 91.2 and 206 GeV 
together with the fitted NNLO QCD predictions. The error bars correspond to statistical uncer- 
tainties. The fit ranges cover the central regions indicated by the solid curves, the theoretical 
predictions extrapolate well outside the fit ranges, as shown by the dotted curves. The plotted 
distributions are scaled by arbitrary factors for presentation. 

NLO+NLLA fits [28]. 

5.2 Theoretical Uncertainties 

In the case of pure fixed-order predictions, the main source of arbitrariness in the predic- 
tions is the choice of the renormalisation scale x^. The residual dependence of the fitted 
value of a s (Mz 2 ) on the renormalisation scale is shown in Fig. 6. A dramatic reduction 
of the scale dependence by a factor of two is observed when going from NLO to NNLO. 
But also compared to NLO+NLLA a net improvement is obtained when fitting the NNLO 
prediction. 

The systematic uncertainty related to missing higher orders is estimated with the 
uncertainty-band method recommended in Ref. [39]. Briefly, this method derives the un- 
certainty of a s from the uncertainty of the theoretical prediction for the event-shape distri- 
bution and proceeds in three steps. First a reference perturbative prediction, here NNLO 
with x^ = 1, is determined using the value of a s obtained from the combination of the 
six variables and eight energies, as explained in section 6. Then variants of the prediction 
with different choices for the x^ scale are calculated with the same value of a s . In each 
bin of the distribution for a given variable, the largest upward and downward differences 
with respect to the reference prediction are taken to define an uncertainty band around 
the reference theory. In the last step, the value of a s in the reference prediction is varied, 
in order to find the range of values which result in predictions lying inside the uncertainty 



- 10 - 



H 
■o 

o 
e 



ALEPH data 
10 r E,. n =91.2 GeV 



10 



fit range 

NNLO, a s =0.1274 ±0.0003, x 2 /n dor =0.44 
NLO, a s =0.1446 ±0.0003, x 2 /n dor =1.7 
NLO + NLLA 

a s =0.1271 ±0.0002, x 2 /n doI =0.62 




0.6 0.65 0.7 0.75 0.8 0.85 0.9 0.95 1 

T 



1 10 

"6 





10 



10 



10 



10 



5 0.2 

iB 

a 



■a 

~-0.2 



-0.4 



ALEPH data 
E =91.2 GeV 



fit range 



NNLO, a =0.1261 ±0.0003, % 7n dor =1.2 
NLO, a,=0.1354 ±0.0003, x 2 'n dd[ =4.3 
NLO + NLLA, a i= 0.1198 ±0.0002, x 2 'n dlll =7.6 



III 
















stat. ® exp. uncertainty 







































statistical uncertainty 




■ , , £ , 








1 



0.05 0.1 0.15 0.2 0.25 0.3 0.35 

M„ 




0.05 0.1 0.15 0.2 0.25 0.3 0.35 0.4 

H.. 




Figure 4: Distributions measured by ALEPH at LEP1, after correction for backgrounds and 
detector effects, of thrust, heavy jet mass, total and wide jet broadening. Fitted QCD predictions 
at different orders of perturbation theory are overlaid. The lower insets show a relative comparison 
of data and QCD fits. 
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C -ln(y 3 ) 

Figure 5: Distributions measured by ALEPH at LEP1, after correction for backgrounds and 
detector effects, of C-parameter and the three-jet transition variable. Fitted QCD predictions at 
different orders of perturbation theory are overlaid. The lower insets show a relative comparison of 
data and QCD fits. 

band for the fit range under consideration. In contrast to the original method [39] we do 
not request the reference prediction to lie strictly inside the uncertainty band, since for the 
present NNLO calculations the latter is still subject to statistical fluctuations. Instead, 
we make a fit of the reference theory with a s as free parameter to the uncertainty band, 
which includes the statistical uncertainty on the C coefficient. The values of a s fitted to 
the upper and lower contour of the uncertainty band finally set the perturbative systematic 
uncertainty. The upward and downward uncertainties are very similar in magnitude and 
the larger is quoted as symmetric uncertainty. The method is illustrated in Fig. 7, taking 
thrust and heavy jet mass as examples. 

The theoretical error depends on the fit range and on the absolute value of a s , scaling 
approximately with at NLO and af at NNLO. This behaviour can be seen from (2.11), 
where variations of x M affect all terms, and individual terms compensate each other across 
the different orders. All compensating terms through a?, but no compensating terms at 
otg are included at NNLO, such that any residual dependence on is to be attributed to 
these terms. 

At LEP2 energies the statistical fluctuations are large. In order to avoid biases from 
downward fluctuations, the theoretical uncertainties are calculated with the value of a s 
obtained by the global combination procedure. For each energy point and in each variable, 



- 12 - 



0.16 



^0.155 using thrust at LEP1 

N 

^ 0.15 - 
0.145 




0.15 



0.145 using heavy jet mass at LEP1 



g" 0.14 



0.135 



0.13 



0.125 



0.12 



0.115 




0.11 



0.5 



1.5 



2 2.5 
x,, 



Figure 6: Dependence of the extracted a s on the renormalisation scale when fitting the thrust 
(left) and heavy jet mass (right) distributions with predictions at different orders of perturbation 
theory. 



the error is evolved to the appropriate energy scale and the uncertainty is calculated for 
the fit range used for the different variables. 

The combined value of a s , used to derive the systematic perturbative error, depends 
itself on the theoretical error. Hence the procedure of calculating the a s combination and 
perturbative error is iterated until convergence is reached, typically after two iterations. 

An additional error is evaluated for the b-quark mass correction procedure. This 
correction has only been calculated to 0{a 2 s )\ no resummed and NNLO expressions are 
yet available. The difference in a s obtained with and without mass corrections is taken as 
systematic error. The total perturbative uncertainty quoted in the tables is the quadratic 
sum of the errors for missing higher orders and for the mass correction procedure. The 
total perturbative error is between 3% and 4% at Mz and decreases to between 2% and 
3% at LEP2 energies. 

The hadronisation model uncertainty is estimated by comparing the standard hadron- 
level event generator programs HERWfG and ARIADNE to PYTHIA for both hadronisa- 
tion and detector corrections. The same set of corrections as in Ref. [28] is used. Both 
corrections are calculated with the same generator in order to obtain a coherent descrip- 
tion at the hadron level. The maximum change with respect to the nominal result using 
PYTHIA is taken as systematic error. It was verified that calculating the detector cor- 
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Figure 7: Theoretical uncertainties for the distributions of thrust (left) and heavy jet mass (right) 
at LEP1. The filled area represents the perturbative uncertainties of the distribution for a given 
value a". The curves show the reference prediction with a® ± Aa s . The theoretical uncertainty 
Aa s is derived from a fit of the reference theory to the contour of the uncertainty band for the 
actual fit range. 



rections and the hadronization corrections with different event generator programs yields 
consistent results. At LEP2 energies the hadronisation model uncertainty is again subject 
to statistical fluctuations. These fluctuations are observed from one energy to the next 
and originate from limited statistics of the fully simulated detector-correction functions. 
Since non-perturbative effects are expected to decrease with 1/Q, the energy evolution of 
hadronisation errors has been fitted to a simple A + B/Q parametrisation. The fit was 
performed for each variable separately. In the fit procedure a weight scaling with luminos- 
ity is assigned to the hadronisation uncertainty at each energy point. This ensures that 
the hadronisation uncertainty at Mz, which is basically free of statistical fluctuations, is 
not altered by the procedure. As in the case of experimental systematic uncertainties, the 
hadronisation uncertainty is essentially identical to that published in [28]. 

The perturbative component of the error, which is the dominant source of uncertainty 
in most cases, is highly correlated between the energy points. The perturbative errors 
decrease with increasing Q, and faster than the coupling constant itself. The overall error 
is in general dominated by the renormalisation scale dependence. 
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6. Combined Results 

The measurements obtained with the six different variables are combined into a single 
measurement per energy using weighted averages. The same procedure as in [28] is applied 
here and it was verified that the combined results using NLO+NLLA with the fit range 
and weights reported in [28] can be reproduced. A weight is assigned to each observable- 
dependent measurement a* proportional to the inverse square of its total error, wi oc l/of. 
The weighted average a s is then given by 



for which the correlation coefficients pij are needed. This correlation between fits of a s to 
different variables is obtained using a large number of simulated data samples and turns 
out to be typically 60%-80% (cf. Ref. [28]). The correlation of systematic errors is taken 
into account by recomputing the weighted average for all variations of the analysis, and 
the change in a s with respect to the nominal value is taken as error. 

The combination of experimental systematic uncertainties at LEP2 energies is obtained 
using a luminosity-weighted average of the uncertainties between 133 GeV and 206 GeV. 
Combined results are given in Table 4, where RMS denotes the root mean squared be- 
tween the different observables. They are shown in Fig. 8, together with a fit of the QCD 
expectation. The curve is seen to be in good agreement with the measurements. In the 
definition of the % 2 of the fit only the uncorrelated component of the errors is taken into 
account, which excludes the perturbative error. 

The combined measurements between 133 and 206 GeV are evaluated at the scale of 
the Z boson mass by using the predicted energy evolution of the coupling constant, eq. 
(2.10). The measurements evolved to Mz are given in Table 5. They are again combined 
using a weighted average, with weights proportional to the inverse square of the total un- 
certainties. In contrast to the combination from different variables, here the measurements 
are statistically uncorrelated. Correlations between systematic uncertainties are taken into 
account and all variations of the determination of a s have been performed for the weighted 
average. 

The final result is a s (Mz 2 ) = 0.1240 ± 0.0033, and the error components are given 
in Table 6. Included in Table 6 is the combination of measurements at LEP2 energies 
without the point at Mz- The total uncertainty of the combined LEP2 measurements is 
smaller than the uncertainty of the LEP1 measurement, because the dominant perturbative 
uncertainties are reduced at higher energies, even after evolution to Mz- 



N=6 




i=l 



and the combined statistical error is 
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Figure 8: The measurements of the strong coupling constant a s between 91.2 and 206 GeV. 
The results using the six different event-shape variables are combined with correlations taken into 
account. The inner error bars exclude the perturbative uncertainty, which is expected to be highly 
correlated between the measurements. The outer error bars indicate the total error. A fit of the 
three-loop evolution formula using the uncorrelated errors is shown. The shaded area corresponds 
to the uncertainty in the fit parameter Aj^jL = 284 ± 14 MeV of the three- loop formula, eq. (2.10). 

7. Discussion 

The improvements achieved by the NNLO fit compared to that at NLO are twofold. First, 
the renormalisation scale uncertainty on the extracted value of a s is significantly reduced. 
This is clearly illustrated in Fig. 6 by comparing the NNLO (solid) and NLO (dashed) 
curves, in the case of thrust and heavy jet mass. Similar improvements are obtained for 
the other variables. Second, we find a better description of the event-shape distributions 
over a larger range. This is evident in Figs. 4 and 5. 

Of course, we also have to perform a careful comparison of the NNLO results with 
those found in NLO+NLLA analyses, which were the state of the art during the LEP era. 
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Figure 9: The measurements of the strong coupling constant a s for the six event shapes, at 
i/s = Mz, when using QCD predictions at different approximations in perturbation theory. 

Once again, Fig. 6 shows that the NNLO perturbative uncertainty is reduced by about 
30% compared to NLO+NLLA. 

It is also remarkable that the a s values obtained from fits to different event shapes 
with NNLO predictions are considerably more self-consistent than those found with either 
NLO or NLO+NLLA expansions. Not only are the extracted values of a s more precise, 
but the spread obtained from the different observables is smaller. This is clearly shown for 
the data set at = Mz in Fig. 9. The key to this dramatic improvement is the rather 
different size of the NNLO corrections to the various observables. 

Despite these improvements our final combined result on a s (M§) still appears to be 
larger than the world average [5]. We recall that the value of a s (M|) obtained from fits 
with NLO+NLLA predictions is smaller than that obtained with pure NLO calculations 
alone. Here we observe that when going from NLO to NNLO there is also a trend in the 
direction of lower values of a s (M|). 

Clearly, resummed predictions are mandatory in the two-jet region. Figures 4 and 
5 clearly show the improvement achieved with NLO+NLLA predictions in the two-jet 
region. Measurements of a s using NLO+NLLA approximations profit from an extended fit 
range in this region. While a consistent matching of NNLLA predictions to NNLO would 
require the analytic resummation of next-to-next-to-leading logarithmic terms, which are 
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not known at present, a matching of existing NLLA expressions to the NNLO calculations 
requires only the calculation of certain matching coefficients [49]. In this case also a better 
description in the two-jet region can be expected. First preliminary results have been 
obtained by us with such matched NNLO+NNLA predictions and fits to data seem to 
confirm the expected trend towards lower values of a s (M§). We will address this issue in 
a forthcoming publication. 

Electroweak corrections may also be of a similar size as the NNLO corrections discussed 
here. At present, these corrections are only known for the inclusive observables like the 
hadronic cross section [50], where they are found to be sizeable [51], but not for the event- 
shape distributions. This issue deserves further study and will also be addressed in a 
forthcoming publication. 

8. Conclusions 

In this paper we used the newly derived NNLO QCD corrections to event shapes in e + e~ 
annihilation [26, 27] to perform the first determination of the strong coupling constant 
a s from event-shape data at NNLO. Our analysis is based on the full set of event-shape 
distributions measured by the ALEPH collaboration [28] at LEP1 and LEP2. 

We observe that the inclusion of NNLO QCD corrections to the different shape vari- 
ables yields several important effects, when compared to the previously available determi- 
nations of a s (M|), based either on pure NLO calculations or NLO predictions matched to 
NLL approximations : 

(a) The dominant theoretical uncertainty on a s (M|), as estimated from scale variations, 
is reduced by a factor 2 (1.3) compared to NLO (NLO+NLLA). A further improve- 
ment can be anticipated from a matching of NNLO and NLLA predictions. 

(b) The central value obtained at NNLO, 

a s (M|) = 0.1240 ± 0.0008 (stat) ± 0.0010 (exp) ± 0.0011 (had) ± 0.0029 (theo), 

is about 10% lower than at NLO, and therefore closer to, albeit still larger than, the 
world average from other observables. It is also larger than the central value obtained 
at NLO+NLLA [28], which shows the obvious need for a matching of NNLO+NLLA 
for a fully reliable result. 

(c) The scatter among the values of a s (M|) extracted from the six different event-shape 
variables is reduced very substantially at NNLO. This is a clear indication that this 
scatter was largely due to missing higher order perturbative corrections in previous 
studies. 

These observations visibly illustrate the improvements gained from the inclusion of 
the NNLO corrections, and highlight the need for further studies on the matching of 
NNLO+NLLA, on the derivation of NNLLA resummation terms, and on the electroweak 
corrections to event shapes in e + e~ annihilation. 
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Table 1: Results for a s (Q) as obtained from NNLO fits to distributions of event-shape variables 
at Q = sfs = 91.2, 133, 161 and 172 GcV. 



Q = 91.2 GeV 


variable 


T 


-lnj/3 


M H 


C 


B w 


B T 


a s 


0.1274 


0.1200 


0.1261 


0.1288 


0.1240 


0.1263 


stat. error 


0.0003 


0.0003 


0.0003 


0.0004 


0.0003 


0.0002 


exp. error 


0.0008 


0.0010 


0.0011 


0.0008 


0.0007 


0.0007 


pert, error 


0.0042 


0.0039 


0.0029 


0.0039 


0.0029 


0.0047 


hadr. error 


0.0019 


0.0016 


0.0044 


0.0017 


0.0018 


0.0018 


total error 


0.0047 


0.0043 


0.0054 


0.0043 


0.0035 


0.0051 


fit range 


0.76-0.90 


1.6-4.4 


0.08-0.18 


0.38-0.72 


0.08-0.15 


0.13-0.28 


Q = 133 GeV 


variable 


T 


-lny 3 


M H 


C 


B w 


Bj' 


a s 


0.1197 


0.1198 


0.1216 


0.1208 


0.1211 


0.1177 


stat. error 


0.0035 


0.0046 


0.0042 


0.0030 


0.0030 


0.0026 


exp. error 


0.0010 


0.0007 


0.0011 


0.0011 


0.0005 


0.0012 


pert, error 


0.0034 


0.0034 


0.0025 


0.0032 


0.0022 


0.0038 


hadr. error 


0.0014 


0.0009 


0.0028 


0.0012 


0.0011 


0.0014 


total error 


0.0052 


0.0058 


0.0058 


0.0047 


0.0039 


0.0049 


fit range 


0.75-0.94 


1.6-4.8 


0.06-0.25 


0.22-0.75 


0.05-0.25 


0.09-0.30 


Q = 161 GeV 


variable 


T 


-lnj/3 


M H 


C 


Bw 


Bx 


a s 


0.1239 


0.1146 


0.1267 


0.1259 


0.1223 


0.1231 


stat. error 


0.0054 


0.0074 


0.0078 


0.0046 


0.0046 


0.0041 


exp. error 


0.0010 


0.0007 


0.0011 


0.0011 


0.0005 


0.0012 


pert, error 


0.0031 


0.0032 


0.0023 


0.0029 


0.0020 


0.0035 


hadr. error 


0.0012 


0.0006 


0.0022 


0.0011 


0.0008 


0.0013 


total error 


0.0064 


0.0081 


0.0085 


0.0056 


0.0051 


0.0056 


fit range 


0.75-0.94 


1.2-5.6 


0.06-0.25 


0.22-0.75 


0.05-0.25 


0.03-0.30 


Q = 172 GeV 


variable 


T 


-lny3 


M H 


C 


Bw 


Bt 


a s 


0.1101 


0.1079 


0.1095 


0.1104 


0.1083 


0.1125 


stat. error 


0.0072 


0.0089 


0.0087 


0.0064 


0.0059 


0.0056 


exp. error 


0.0010 


0.0009 


0.0011 


0.0011 


0.0005 


0.0012 


pert, error 


0.0030 


0.0032 


0.0022 


0.0028 


0.0019 


0.0033 


hadr. error 


0.0012 


0.0005 


0.0021 


0.0010 


0.0007 


0.0012 


total error 


0.0079 


0.0093 


0.0092 


0.0069 


0.0063 


0.0067 


fit range 


0.80-0.96 


1.2-5.6 


0.06-0.25 


0.22-0.75 


0.05-0.25 


0.09-0.30 
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Table 2: Results for a s (Q) as obtained from NNLO fits to distributions of event-shape variables 
at g = V« =183, 189, 200 and 206 GeV. 



Q = 183 GeV 


variable 


T 


— lm/3 


M H 


C 


B w 


B T 


a s 


0.1132 


0.1070 


0.1133 


0.1160 


0.1114 


0.1117 


stat. error 


0.0032 


0.0041 


0.0039 


0.0029 


0.0027 


0.0027 


exp. error 


0.0011 


0.0009 


0.0012 


0.0014 


0.0005 


0.0012 


pert, error 


0.0029 


0.0031 


0.0022 


0.0027 


0.0019 


0.0033 


hadr. error 


0.0011 


0.0004 


0.0019 


0.0010 


0.0007 


0.0012 


total error 


0.0046 


0.0052 


0.050 


0.0043 


0.0033 


0.0046 


fit range 


0.75-0.94 


1.2-5.6 


0.06-0.25 


0.22-0.75 


0.05-0.25 


0.09-0.30 


Q = 189 GeV 


variable 


T 


-lny 3 


M H 


C 


B w 






0.1140 


0.1077 


0.1118 


0.1122 


0.1093 


0.1124 


stat. error 


0.0020 


0.0025 


0.0025 


0.0017 


0.0016 


0.0015 


exp. error 


0.0010 


0.0007 


0.0013 


0.0011 


0.0005 


0.0013 


pert, error 


0.0028 


0.0031 


0.0021 


0.0026 


0.0018 


0.0030 


hadr. error 


0.0011 


0.0004 


0.0018 


0.0009 


0.0006 


0.0012 


total error 


0.0037 


0.0041 


0.0040 


0.0034 


0.0025 


0.0038 


fit range 


0.80-0.96 


1.6-5.6 


0.04-0.20 


0.18-0.60 


0.04-0.20 


0.075-0.25 


Q = 200 GeV 


variable 


T 




M H 


C 


Bw 


Bt 


Ot s 


0.1094 


0.1071 


0.1079 


0.1105 


0.1071 


0.1082 


stat. error 


0.0022 


0.0026 


0.0027 


0.0018 


0.0017 


0.0017 


exp. error 


0.0010 


0.0007 


0.0014 


0.0011 


0.0005 


0.0013 


pert, error 


0.0027 


0.0031 


0.0021 


0.0026 


0.0018 


0.0029 


hadr. error 


0.0010 


0.0003 


0.0017 


0.0009 


0.0006 


0.0011 


total error 


0.0037 


0.0041 


0.0041 


0.0035 


0.0026 


0.0038 


fit range 


0.80-0.96 


1.6-5.6 


0.04-0.20 


0.18-0.60 


0.04-0.20 


0.075-0.25 


Q = 206 GeV 


variable 


T 


-lnt/3 


M H 


C 


Bw 


Bj 1 


a s 


0.1075 


0.1023 


0.1093 


0.1064 


0.1066 


0.1056 


stat. error 


0.0021 


0.0026 


0.0025 


0.0017 


0.0016 


0.0016 


exp. error 


0.0010 


0.0007 


0.0011 


0.0011 


0.0005 


0.0012 


pert, error 


0.0026 


0.0031 


0.0021 


0.0025 


0.0017 


0.0029 


hadr. error 


0.0010 


0.0003 


0.0016 


0.0009 


0.0005 


0.0011 


total error 


0.0037 


0.0041 


0.0038 


0.0034 


0.0025 


0.0037 


fit range 


0.80-0.96 


1.6-5.6 


0.04-0.20 


0.18-0.60 


0.04-0.20 


0.075-0.25 
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Table 3: Various fit results at a centre-of-mass energy of 91.2 GeV, using different predictions of 
perturbative QCD, with the renormalisation scale fixed to fj, = Mz- 





T 


C 


M H 


B w 


Bj' 


-lny 3 


fit range 


0.76 - 0.90 


0.38 - 0.72 


0.08 - 0.18 


0.08 - 0.15 


0.13 - 0.28 


1.6-4.4 


NNLO 


0.1274 


0.1288 


0.1261 


0.1240 


0.1263 


0.1200 


X 2 /N dof 


5.73/13 


9.07/16 


10.8/9 


3.89/6 


22.6/14 


14.5/13 




= 0.44 


= 0.57 


= 1.20 


= 0.65 


= 1.61 


= 1.11 


stat. error 


0.0003 


0.0004 


0.0003 


0.0003 


0.0002 


0.0003 


NLO 


0.1446 


0.1437 


0.1353 


0.1323 


0.1454 


0.1246 


X 2 /N dof 


1.72 


1.4 


4.3 


1.7 


11.9 


2.1 


NLO+NLLA 


0.1271 


0.1235 


0.1198 


0.1164 


0.1269 


0.1180 


X 2 /N dof 


0.62 


1.66 


7.6 


7.61 


2.2 


5.4 



Table 4: Combined results for a s (Q) using NNLO predictions. 



Q[GeV] 


91.2 


133 


161 


172 


183 


189 


200 


206 


a s (Q) 


0.1252 


0.1202 


0.1231 


0.1099 


0.1122 


0.1110 


0.1082 


0.1064 


stat. error 


0.0002 


0.0023 


0.0036 


0.0048 


0.0021 


0.0013 


0.0014 


0.0013 


exp. error 


0.0008 


0.0009 


0.0009 


0.0009 


0.0009 


0.0009 


0.0009 


0.0009 


pert, error 


0.0037 


0.0030 


0.0028 


0.0027 


0.0026 


0.0024 


0.0024 


0.0023 


hadr. error 


0.0018 


0.0012 


0.0009 


0.0009 


0.0008 


0.0008 


0.0007 


0.0007 


total error 


0.0042 


0.0041 


0.0047 


0.0057 


0.0035 


0.0030 


0.0030 


0.0029 


RMS 


0.0031 


0.0014 


0.0043 


0.0017 


0.0030 


0.0023 


0.0013 


0.0023 



Table 5: Combined results for a s (Mz) using NNLO predictions. 



Q[GeV] 


91.2 


133 


161 


172 


183 


189 


200 


206 


a s {M z ) 


0.1252 


0.1276 


0.1352 


0.1207 


0.1247 


0.1238 


0.1215 


0.1196 


stat. error 


0.0002 


0.0026 


0.0043 


0.0058 


0.0026 


0.0016 


0.0017 


0.0017 


exp. error 


0.0008 


0.0011 


0.0012 


0.0012 


0.0012 


0.0011 


0.0011 


0.0011 


pert, error 


0.0037 


0.0032 


0.0030 


0.0030 


0.0028 


0.0027 


0.0027 


0.0027 


hadr. error 


0.0018 


0.0013 


0.0011 


0.0011 


0.0010 


0.0009 


0.0009 


0.0009 


total error 


0.0042 


0.0045 


0.0055 


0.0067 


0.0041 


0.0035 


0.0035 


0.0035 


RMS 


0.0031 


0.0015 


0.0053 


0.0020 


0.0037 


0.0029 


0.0017 


0.0029 
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Table 6: Weighted average of combined measurements for a s (Mz) obtained at energies from 91.2 
GeV to 206 GeV and the average without the point at y/s = Mz ■ 



data set 


LEP1 + LEP2 


LEP2 


a s (M z ) 


0.1240 


0.1238 


stat. error 


0.0008 


0.0009 


exp. error 


0.0010 


0.0011 


pert, error 


0.0029 


0.0028 


hadr. error 


0.0011 


0.0010 


total error 


0.0033 


0.0033 
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